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We solve the time-dependent Gross-Pitaevskii in 2-D 
to simulate the flow of an object through a dilute Bose 
condensate trapped in a harmonic well. We demon- 
strate vortex formation and analyze the process in 
terms of the accumulation of phase-slip and the evolu- 
tion of the fluid velocity. 
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The origin of drag in quantum liquids is central 
to the understanding of superfluidity Q. For a 
weakly-interacting fluid, Bogoliubov showed that 
macroscopic occupation of the ground state leads 
to a linear dispersion curve, and hence superfluid- 
ity for motion slower than a critical velocity. How- 
ever, for liquid helium, the transition to normal 
flow is observed at a much lower velocity than can 
be explained by the dispersion curve. This led 
Feynman |^ to suggest that the onset of dissipa- 
tion may arise due to vortex shedding, but exper- 
imental verification of this idea has been impeded 
by the two-fluid nature of liquid helium, which 
complicates quantitative comparison with theory. 

The recent experimental discovery of Bose conden- 
sation in dilute alkali vapours presents a near- 
perfect system for the study of quantum fluids: 
the condensates are almost pure, and sufficiently 
dilute that the interactions can be accurately pa- 
rameterized in terms of a scattering length. As a 
result, a relatively simple non-linear Schrodinger 
equation (NLSE), known as the Gross-Pitaevskii 
equation Ig] , gives a precise description of conden- 
sate dynamics. Experiments have confirmed that 
the NLSE is remarkably accurate in the limit of 
low temperature [R-|9| . 



Recent experiments have also demonstrated that 
far-off resonant laser light may be used to split 
pO| or excite the condensate The question 

arises whether such light forces may be used to 
study the phenomenon of vortex formation, and 
thereby shed some light on the issue of superflu- 
idity. The appearance of vortices in solutions of 
the 2-D NLSE has been considered by Frisch et al. 
Jia. However, they consider homogeneous fluid 



flow past an impenetrable obstacle, which is inap- 
propriate for atomic condensates. Here, we solve 
the time-dependent 2-D NLSE to simulate the flow 
of an 'object' through a dilute inhomogeneous con- 
densate. The 'object' is the potential barrier pro- 
duced by far-off resonant blue-detuned laser beam 
P|JlO[|. The condensate is formed in a harmonic 
trap with the obstacle at the center. Subsequently, 
the obstacle is withdrawn at a constant speed, v 
(from the Galilean invariance of the NLSE W,this 
is equivalent to ffow past a stationary objectjl Vor- 
tices are formed when the condensate, whose mo- 
tion is limited by the speed of sound, cannot adjust 
sufficiently quickly to the density wave induced by 
the object. We show how the time scale for vortex 
formation can be extracted from the phase-slip of 
the wavefunction. 

The evolution of the condensate is described by the 
NLSE [ITl: 



idtip^ (-v^ + y + cii/'P) V' 



(1) 



where C — 8nNa is the nonlinear coefficient, N is 
the number of atoms in the condensate, and a is 
the s-wave scattering length. The potential term, 
V = j{x^ + y2) + aexp{-f3x^ - (3{y - vtf), de- 
scribes a symmetric harmonic trap with a light- 
induced 'Gaussian' potential barrier, moving with 
velocity, v. We have considered a variety of ob- 
stacle parameters but present detailed results for 
C = 500, a = 30, and /3 = 3. 

The eigenvalue problem, ?/; = ^(a;, j/)e~'''*, (where 
/i is the chemical potential), and w = 0, is easily 
solved using finite difference methods. With a = 0, 
we find ^ = 9.003. The depletion of density in the 
trap center due to the presence of the obstacle, 
(with a = 30, /? = 3, w = 0), raises the chemical 
potential to = 9.208. A change in the object 
potential creates a local disturbance which prop- 
agates through the fluid. The time evolution of 
"0 was determined by the split-operator technique 



14 . Extinction of the object, excites sound waves 
11 , which propagate with speed, c = ^J2C\%l)\'^, 
(^ 4.3 at the peak density). Translation of the ob- 
ject displaces the condensate center-of-mass p^ , 
and leads to vortex-pair creation as illustrated in 
Fig. 1^. The vortex minima are highlighted by the 
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surface plot shown in Fig. ||. The half-width of the 
vortex core is comparable to the healing length, 
A ~ C~^l'^ . We have simulated the free expansion 
of the condensate and observe that the vortices also 
expand, permitting experimental detection using 
optical imaging. 
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FIG. 1. Contour plot of the condensate density, 
[■(/'(a;, t/, 3.0)p for v — 2.0. There are 25 contours, 
equally spaced between 0.0 and 0.0182. 




FIG. 2. Surface plot of the condensate density, 
J/, 3.0)p for v = 2.0, viewed from below, illus- 
trating the minima corresponding to the object and 
the vortex-pair. 

For low object speeds, no vortices are observed 
implying that the process requires motion faster 
than a critical speed. For an impenetrable cylin- 
der of diameter d, the critical speed is related 



to the energy required to create a vortex-pair, 
e — (C/2a) ln(d/A) For a finite potential bar- 
rier, we find that the vortices emerge from a point, 
i.e., initially d '-^ A, so the pair creation energy, 
e ^ 0, and the critical speed depends principally on 
the speed of sound. However, due to the conden- 
sate inhomogeneity, the speed of sound and hence 
the critical speed are spatially dependent. For this 
reason, we choose to study the process of vortex 
formation in terms of the fluid velocity and the 
phase of the wavefunction. 

Fig. H shows a quiver plot of the fluid velocity, 
■*^s = (■0*^^ ~ in the vicinity of the 

object, for t — 3.0 and v = 2.0 (as in Figs. |l| and 
One sees that vortices are formed in pairs with op- 
posing vorticity, and that the velocity is inversely 
proportional to distance from the vortex line, as 
expected. Two pairs are discernable in Fig. ^ one 
&i y = 3.5, which corresponds to the density zeros 
in Fig. 0, and a second at y = 5.9, (not visible 
in Fig. PI because it has yet to separate from the 
object). The second pair appears due to the rapid 
accumulation of phase-slip at the edge of the con- 
densate, as discussed below. 
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FIG. 3. Plot showing the fluid velocity, Vs, in the 
vicinity of the object, for t = 3.0 and v = 2.0. The 
arrow length is 40ws in h.o.u. 

By calculating the wavefunction phase, S{x,y,t), 
one may verify that the phase changes by 27r on 
circulation of the vortex line In addition, the 
evolution of the phase determines the time-scale for 
vortex formation. At t = 0.0, the phase is uniform, 
then the motion induces a dephasing or 'phase-slip' 
centred on the object. The phase-slip between two 
neighbouring points along the y-axis is defined as: 
ASiy, t) = arg {^'(O, y + Ay,t)}- arg {^'(O, y, t)}, 
with — TT < arg {z} < tt. The maximum phase-slip, 
AS'max, occurs at a y-coordinate close to the center 
of the object. 
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Fig. H shows a plot AS'^ax for different object 
speeds. For v — 2.0 and 3.0, the phase-slip ac- 
cumulates gradually, reaching a value of tt, and 
then changing sign at a critical time, tc- The sign 
change indicates a reversal of the flow, and may be 
taken to define the moment of vortex creation. The 
critical time is proportional to the speed of sound: 
vortices are produced faster for smaller c, as the 
condensate requires longer to adjust to the pertur- 
bation (compare curves for C ~ 200 and C — 500 
in Fig. For t > tc, the phase-slip builds again 
until another vortex-pair is formed. For v = 1.3, 
(upper plot in Fig. ^) the phase-slip initially sat- 
urates and then increases towards the edge of the 
condensate, where the density and hence the local 
speed of sound are lower. This explains the sud- 
den formation of a vortex-pair as the object leaves 
the condensate, (e.g. the second pair for v — 2.0 
appearing in Figs. g| and ^). For v = 3.0, a third 
pair is created a,t t — 2.5: the object is still sur- 
rounded by condensate at this time because of a 
motion induced 'stretching'. 





direction of motion, as the fluid attempts to fill the 
'hole' left by the departing object. At i = 1.1, the 
fluid begins to skirt around the object and fill the 
hole with a forward flow. The wavefunction nodes 
(separated by less than one grid point at this stage) 
are pulled apart due to the gradient in the object 
potential. By t — 1.2, their separation is compara- 
ble to the healing length, and the pattern of vortex 
flow is beginning to emerge. For t > 1.2 the vortex 
pair continue to separate, and the on-axis forward 
flow gradually decreases, reverting to a reverse flow 
at t = 2.0 (i.e., when the phase slip changes sign, 
see Fig. Thereafter, the process repeats. For 
longer times, the motion of the vortex pair is sen- 
sitive to the exact parameters used, and a detailed 
study is required before any generalizations can be 
made. 
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FIG. 5. Plot showing the velocity field in the vicinity 
of the object for v — 2.0, at times spanning the instant 
of vortex formation: t — 1.0 (top), 1.1 (middle), and 
1.2 (bottom). The arrow length is lOvs in h.o.u. 
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FIG. 4. Time-evolution of the maximum phase-slip, 
AS'max, for V = 1.3 (top); V = 2.0 (middle); and v — 3.0 
(bottom). For v = 2.0 we show an additional curve for 
C = 200 (dotted), illustrating that vortices are pro- 
duced faster when the speed of sound is slower. 

The evolution of the velocity field while the vortex 
is forming is shown in Fig. |5|. For t = 1.0 {t < tc), 
the on-axis flow is backwards, i.e., opposite to the 



To summarize, we have solved the time-dependent 
NLSE in 2-D to simulate the flow of an object 
through a dilute Bose condensate trapped in a 
harmonic well. The process of vortex formation 
has been analyzed in terms of the accumulation of 
phase-slip and the evolution of the velocity field. 
We find that the vortices emerge from a point close 
to the center of the object, and that the vortex 
shedding frequency is higher when the speed of 
sound is lower. 
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A complete deseription of dilute atomic vapours 
should include the effects of the non-condensate 
density, finite temperature, and dissipation. To 
what extent the NLSE provides an accurate de- 
scription of vortex formation requires quantitative 
comparison between simulation and experiment. 
Such a comparison will provide the focus for fu- 
ture work. 
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m, the units of length, time and energy are 
{h/2mLu)^^^ , Lu^^ and huj, respectively. 

[14] Most of the calculations were performed using a 
square box of side 10 h.o.u. divided into a grid of 
256x256 points. A larger grid (20 h.o.u., 512x512) 
was used to check for edge-effects. More details 
of the numerical methods employed arc discussed 
elsewhere: B. Jackson, J. F. McCann, and C. S. 
Adams, to be published. 

[15] The center-of-mass motion provides a useful di- 
agnostic: the initial acceleration is v divided by 
the integration time-step; after a few steps the so- 
lution is only sensitive to v; and when the object 
leaves the condensate, the center-of-mass oscillates 
at the trap frequency. 
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